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In this paper we present the results of a coherent narrow-band search for continuous gravitational- 
wave signals from the Crab and Vela pulsars conducted on Virgo VSR4 data. In order to take into 
account a possible small mismatch between the gravitational wave frequency and two times the 
star rotation frequency, inferred from measurement of the electromagnetic pulse rate, a range of 
0.02 Hz around two times the star rotational frequency has been searched for both the pulsars. No 
evidence for a signal has been found and 95% confidence level upper limits have been computed 
both assuming polarization parameters are completely unknown and that they are known with some 
uncertainty, as derived from X-ray observations of the pulsar wind torii. For Vela the upper limits 
are comparable to the spin-down limit, computed assuming that all the observed spin-down is due 
to the emission of gravitational waves. For Crab the upper limits are about a factor of two below 
the spin-down limit, and represent a significant improvement with respect to past analysis. This is 
the first time the spin-down limit is significantly overcome in a narrow-band search. 


I. INTRODUCTION 

Continuous gravitational- wave signals (CW) are emit- 
ted by spinning neutron stars if asymmetric with respect 
to the rotation axis. The asymmetry can be due to var- 
ious mechanisms, like a non-axismmetric residual strain 
from the star’s birth or a strong internal magnetic field 
not aligned to the rotation axis, see e.g. [l|. When the 
source parameters, position, frequency, spin-down, are 
known with high accuracy targeted searches can be done 
using optimal analysis methods, based on matched filter- 
ing. This happens, for instance, with known pulsars for 
which accurate ephemerides are often obtained from elec- 
tromagnetic (EM) observations, especially in the radio, 
gamma-ray and X-ray band. This means that a strict cor- 
relation between the gravitational wave frequency /o and 
the star’s measured rotational frequency f rot is assumed. 
In the classical case of a non-axisymmetric neutron star 


rotating around one of its principal axes of inertia the 
gravitational frequency would be exactly twice the rota- 
tion frequency of the star. Several targeted searches for 
CW have been conducted on data from first generation 
interferometric detectors. No evidence for a signal has 
been found, but interesting upper limits have been placed 
in a few cases M. Given the uncertainties on gravita- 
tional wave emission mechanisms and also the lack of a 
full detailed picture of the electro-magnetic emission ge- 
ometry, it is not obvious at all that the gravitational- wave 
emission takes place at exactly twice the star measured 
pulse rate, or that such relation holds for observation 
times of months to years. For instance, if a neutron star 
is made of a crust and a core rotating at slightly different 
rates, and if the gravitational- wave emission is dominated 
by an asymmetry in the core then a search targeted at 
2 f rot would assume a wrong signal frequency. We then 
consider here that the signal frequency can be slightly 
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different with respect to fo = 2 f rot and can vary in the 
interval 


f(t ) G [/o(t)(l - 5), f 0 (t)(l + <5)] (1) 

where S is a small positive real number. Following the 
discussion in [5], if the star crust and core are linked 
by some torque that tends to enforce corotation on a 
timescale r c , then S ~ t c / t sc [ , where r s d ^ fo/fo is the 
characteristic spin-down time. A relation of the form of 
Eq. (pQ) also holds in the case the gravitational radiation 
is produced by free precession of a nearly bi-axial star |[(| , 
in which case S is of the order of ( I zz — I ra )// ra , where 
I xx and I zz are the star moments of inertia with respect 
to a principal axis on the equatorial plane and aligned 
with the rotation axis respectively. In general, a value of 
S of the order of, say, 10 -4 , corresponds to r c ^ 10~ 4 t sc i 
which, depending on the specific targeted pulsar can be 
several months or years. This would be comparable or 
larger than the longest timescale observed in pulsar glitch 
recovery where a recoupling between the two components 
might occur. In terms of free precession, 5 ~ 10 -4 is on 
the high end of the range of deformations that neutron 
stars could sustain M- 

Narrow-band searches have not received much atten- 
tion until now, one notable exception being the Crab 
narrow-band search over the first 9 months of LIGO S5 
data Q, based on the ^-statistic 0- In previous work 
m an optimal method, based on matched filtering in the 
space of signal Fourier components, has been proposed 
to carry out a search for CW over a small frequency and 
spin-down range around the values inferred from EM ob- 
servations. In this paper we describe the application of 
such a method to a narrow-band search of CW from Crab 
and Vela pulsars in the data of Virgo VSR4 run. As no 
evidence for a signal is found, we place upper limits on 
signal amplitude. 

The outline of the paper is the following. In Sec.HTlthe 
expected gravitational wave signal from spinning neutron 
stars is introduced. In Sec. umt, he main characteristics of 
Virgo VSR4 data are presented. In Sec. [IV] the analysis 
method is described, while in Sec. |V] the analysis results 
are discussed. Sec. lVIl is dedicated to the validation tests 
of the analysis procedure. Finally, Sec. IVIII contains the 
conclusions. 


II. THE GRAVITATIONAL WAVE SIGNAL 

A non-axisymmetric neutron star steadily spinning 
about one of its principal axis emits a quadrupolar 
gravitational- wave signal at twice the star rotational fre- 
quency f rot , which at the detector can be described [12] 
as 


h(t) = H 0 (H+A+(t) + H X A X (t)yMt)t+<s> 0 ) ( 2 ) 

where taking the real part is understood. The constant 
4>o is the initial signal phase. The angular signal fre- 


quency uo(t) = 47r/ ro t(£), is a function of time. Conse- 
quently the signal phase is not that of a simple monochro- 
matic signal and depends on both the intrinsic rotational 
frequency and frequency derivatives of the pulsar and on 
Doppler and propagation effects, which include the Ein- 
stein delay and, possibly, the Shapiro delay. These varia- 
tions are corrected in the time-domain in a way described 
in Sec. UYl The two complex amplitudes H+ and H x are 
given respectively by 


H+ = 


cos 2'i/j — irj sin 2^ 
+ Trp 


( 3 ) 


sin + vr\ cos 2^ 

yi + 'if 


( 4 ) 


in which r] is the ratio of the polarization ellipse semi- 
minor to semi-major axis and the polarization angle if 
defines the direction of the major axis with respect to the 
celestial parallel of the source (measured counterclock- 
wise). The functions A + and A x describe the detector 
response as a function of time and are a linear combi- 
nation of terms depending on the Earth sidereal angular 
frequency, Q 0 : 


A + = ao + aic cos Q 0 £ + ais sin D 0 £+ 


d 2 c cos2Q 0 £ + d 2 s sin2D 0 £ (5) 


A x = b\ c cos D 0 £ + bi s sin D 0 £+ 

& 2 C cos2D 0 t + b 2 S sin2D 0 £ (6) 

where the coefficients depend on the source position and 
detector position and orientation on the Earth 0 . 

As discussed e.g. in Q, the overall wave amplitude Hq 
and r] are related to the “standard” signal amplitude 

47 T 2 G I zz £f$ 

'*> = ^4 — P) 

and to the angle i between the star rotation axis and the 
line of sight to the source by 


H 0 — ho 


1 + 6 cos 2 l + cos 4 l 


(8) 


2 cos l 
1 + cos 2 i ’ 


( 9 ) 


In Eq. 0 G is the gravitational constant, c is the light 
velocity, I zz is the star moment of inertia with respect 
to the principal axis aligned with the rotation axis, e = 
ixx-iyy ^ e equatorial ellipticity expressed in terms of 

principal moments of inertia, and d is the source distance. 
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Equating the gravitational-wave luminosity to the ki- 
netic energy lost as the pulsar spins down gives us the 
so-called spin- down limit on gravitational- wave strain 


h 


sd 

0 


/ 5 GIzzfvot \ 1/2 
y 2 c 3 d 2 f rot J 


7-1/2 

8.06 x 1CT 19 -^- 

^kpc 


|/rotiy /2 //rot\~ 1/2 

Hz/s J V Hz / 


(10) 


The low-frequency sensitivity of VSR4 was significantly 
better than that of previous Virgo runs, especially due to 
the use of monolithic mirror suspensions, and basically 
in agreement with the planned sensitivity of the initial 
Virgo interferometer. In Fig. [T| a typical VSR4 strain 
sensitivity curve is shown. In Fig. [2] the average power 
spectrum around Crab and Vela reference frequency is 
plotted. Note the large sensitivity improvement around 
the Vela frequency after removal of an instrumental dis- 
turbance, about one month after the beginning of the 
run, see the figure caption for more details. 


where Iss is the star’s moment of inertia in units of 
10 38 kgm 2 , / rot is the time derivative of the rotational 
frequency and dk pc is the distance to the pulsar in kilo- 
parsecs. The spin-down limit on the signal amplitude 
corresponds to an upper limit on the star’s fiducial ellip- 
ticity 


e sd = 0.237 


fr s o d 

10 -24 


f~ 2 T - 1 
J rot i 38 


dki 


pc- 


(ii) 


This quantity, for a given neutron star equation of state, 
can be related to the physical ellipticity of the star sur- 
face [8]. Setting a gravitational-wave upper limit below 
the spin-down limit is an important milestone in CW 
searches as it allows us to constrain the fraction of spin- 
down energy due to the emission of gravitational-waves 
, which gives insight into the spin-down energy budget 
of the pulsar. On the other hand the l = m = 2 mass 
quadrupole moment Q 22 is related to the gravitational- 
wave amplitude by 


Q 22 


, ( 

= 1 / — 1 ^£ = h 0 I 



\167r 2 G/ I 


see e.g. Q, and is independent of any assumptions about 
the star’s equation of state and its moment of inertia. 


III. INSTRUMENTAL PERFORMANCE IN 
VIRGO VSR4 RUN 

Interferometric gravitational wave detectors, LIGO 
0 , Virgo 0 , GEO 0 , have collected a large amount 
of data in recent years. For the analysis described in this 
paper we have used calibrated data, sampled at 4096 Hz, 
from Virgo VSR4 science run. The run extended from 
June 3rd, 2011 (10:27 UTC) to September, 5th 2011 
(13:26 UTC), with a duty factor of about 81%, corre- 
sponding to an effective duration of 76 days. Calibra- 
tion uncertainties amounted to 7.5% in amplitude and 
(40 + 50 fkHz) mrad in phase up to 500 Hz, where fkHz is 
the frequency in kilohertz. The uncertainty on the ampli- 
tude will contribute to the uncertainty on the upper limit 
on signal amplitude, together with that coming from the 
finite size of the Monte Carlo simulation used to compute 
it, see Sec. m The calibration error on the phase can be 
shown to have a negligible impact on the analysis @- 


IV. SEARCH DESCRIPTION 

The analysis pipeline consists of several steps, schemat- 
ically depicted in Fig.[3j which we summarize here. More 
details are given in [11]. The starting point is a collection 
of windowed and interlaced (by half) “short” (1024 sec- 
onds) FFTs (the Short FFT Database - SFDB) built from 
calibrated detector data [ll|. At this stage a first clean- 
ing procedure is applied to the data in order to remove 
big and short time duration disturbances, that could re- 
duce the search sensitivity. A small frequency band is 
then extracted from the SFDB, but large enough to con- 
tain the Doppler and spin-down variations of the signal 
possibly emitted by the target pulsar. In the analyses 
described in this paper, for instance, it was of 0.15 Hz. 
At this point we go back to the time domain (the sam- 
pling frequency is still the original one, 4096 Hz) and 
make barycentric and spin-down corrections. Due to the 
Doppler effect the received frequency f(t) is related to 
the emitted frequency fo(t) by the well-known relation 
(valid in the non-relativistic approximation) 

< 13 > 

where 4>(t) is the signal phase, v = Vrev+Vrot is the detec- 
tor velocity with respect to the Solar system barycenter 
(SSB), given by the sum of the Earth’s orbital velocity, 
v rev , and rotational velocity, v ro u while h is the versor 
identifying the source position. In practice the Doppler 
effect is efficiently corrected in the time domain by chang- 
ing the time stamp t of data samples according to 

Tl = t+ r St± =t + Ar (14) 

c 

where r is the detector position in the SSB. The correc- 
tion term A#, called Roemer delay , amounts up to about 
500 seconds and corresponds to the time taken by a signal 
traveling at the speed of light to cover the distance be- 
tween the detector and the SSB. Unlike Eq. (fl3| ) . Eq. (fl4|) 
does not depend on the frequency, which means that one 
single correction holds for every frequency. In fact there 
are other smaller relativistic effects that must be taken 
into account when making barycentric corrections. One 
is the Einstein delay , A#, due to the Earth motion and 
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FIG. 1. A typical sensitivity curve of Virgo VSR4 run, expressed in terms of noise spectral density. 




FIG. 2. Average spectrum of VSR4 data around Vela (top plot) and Crab (bottom plot) reference frequency (identified by the 
vertical black lines). For Vela two spectra are in fact shown in the same plot. The upper one (black curve) has been computed 
over the first 29 days of the run when the frequency region around Vela was affected by an instrumental disturbance. This 
was due to a non-linear coupling between the DARM control line at 379 Hz and calibration lines at 356 Hz and 356.5 Hz and 
was moved away from the Vela region by shifting of 5 Hertz the frequency of calibration lines. The lower plot (blue curve) is 
the average spectrum computed after the removal of the disturbance. The large spectral disturbance appearing near the Crab 
reference frequency is due to a poorly understood coupling between a line at 60 Hz, part of a 10-Hertz comb of lines of likely 
magnetic origin, and the fundamental pendulum frequency of the Virgo mirror system, at 0.594 Hz. 
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FIG. 3. Scheme of the narrow-band analysis pipeline. The starting point is constituted by detector calibrated data, sampled at 
4096 Hz. After bary centric and spin-down corrections the data are down-sampled at 1 Hz. The number of different spin-down 
corrections applied to the data, ny, and then the number of produced corrected time series, is given by Eq. m- See text and 
0 for more details. 


the gravitational redshift at the Earth geocenter due to 
the solar system bodies, and that amounts to about 2 
milliseconds at most. Another effect is the Shapiro delay 
A s which takes into account the deflection of a signal 
passing near a massive body, and which can be shown 
to be negligible for CW searches, unless the source line 
of sight passes very near the Sun’s limb. See e.g. [0 
for explicit expressions of A# and A s- Overall, we can 
make the full barycentric corrections by introducing the 
new time variable 

n = t + A r + A e ~ A s (15) 

This transformation corresponds to referring the data 
collected at the detector site to the SSB, which can be 
considered an inertial reference frame to a very good ap- 
proximation. In practice, for the Crab and Vela anal- 
yses described in this paper the Shapiro delay can be 
neglected. 

For given values of the signal frequency / and fre- 
quency derivatives /, / we could take into account the 
spin-down in a similar way, by re-scaling time according 
to 


rapidly spinning pulsars demonstrate deviations from a 
Taylor expansion when spinning down, a phenomenon 
known as timing noise. Given our lack of knowledge of 
the exact mechanism that might cause the gravitational 
wave frequency to deviate from (twice) the observed elec- 
tromagnetic frequency, we cannot be sure if we should 
expect the timing noise to be present in the signal we 
are searching for (see (0 for discussion). However, a 
study based on the monthly Crab ephemeris data [18j 
has shown that in the Crab, one of the noisiest of pul- 
sars, the effect of timing noise over the duration of the 
observation period is negligible, producing a mismatch 
between a Taylor expansion and a signal based on the 
actual ‘noisy’ time series of less than 1%. This is small, 
confirming that timing noise is likely to have a negligible 
effect on our analysis. 

In practice the spin-down correction is applied after 
barycentric corrections, then the time t that appears in 
Eq. (fl6l) is in fact the re-scaled time r\ of Eq. (p~5|) . In 
this case the correction depends explicitly on the search 
frequency /. The number of frequency bins which cover 
the range A / = 2/o5, corresponding to Eq. Q, is 


T 2 =t + (16) 

where to is the initial time of the data set and higher 
order terms can be included if needed. 

Note that re-scaling the time in this way to make the 
signal monochromatic assumes that the intrinsic phase 
evolution of the pulsar’s GW signal is described by a Tay- 
lor expansion over the entire observation period. Elec- 
tromagnetic observations of pulsars show that young, 


-VJ) (JL) 

0.02 Hz y \ 1 yr y 
(17) 

with 5f = being the frequency spacing, with T the 
total observation time, and where [ • ] stands for the 
nearest integer. Similarly, we take the width of the first 
order spin-down range as A / = 2|/o|5, and A / = 2\fo\S 
for the second order. Let us consider for the moment 
only the first order term. For a fixed value of / we do 
not want to make an explicit spin-down correction for 


n f = 


AT 

jf. 


= [A/ • T] m 6.3 x 10 s 
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each value of /. So, we fix the value of the frequency in 
the denominator of Eq. (®D, for instance at the reference 
frequency /o, and exploit the fact that the same correc- 
tion that holds for the pair (/o, /o) is also valid for all the 
pairs (/, /) such that 


fo _ / 
/o / 


(18) 


In practice, this means that at each frequency / we ex- 
plore a range of first order spin-down values [/ — , / + 

%£], with 


/-./ 07 - (19) 

Jo 

The number of bins in the spin-down term of the first 
order can be computed by considering the “natural” step 
Sf = 2^2 and is given by 


«400 


n (1) = nj = [2 A f-T 2 

k ) (—) (JL ) 2 

10 10 Hz/s J \10 3 J \1 yr / 


( 20 ) 


Following the same reasoning for the second order spin- 
down we arrive at the following expression for the number 
of steps that should be taken into account: 


~ 0.6 


n ^ m rij = |^A / • T 3 

k ) (JL ) 3 

10 20 Hz/s 2 J \10 3 ) \1 yr / 


(21) 


For observation times of the order of a year, reasonable 
values of second order spin-down, and values of 5 typical 
of narrow-band searches, nj < 1 , which means that we do 
not need to consider more values other than the “central” 
one. In Tab. [T] the Crab and Vela EM-inferred positional 
and rotational parameters are shown. Also estimations 
of the polarization parameters are given, which are used 
in the computation of upper limits, see Sec. [Vl In Fig. 0] 
the portion of the / — / plane actually covered in the 
Crab and Vela narrow-band search is shown. 

Once barycentric and spin-down corrections have been 
done the data are down-sampled to a much lower rate 
with respect to the original one, 1 Hz in the present case. 
This strongly reduces the amount of data to be handled 
in the next step of the analysis. At this stage for a given 
source we have rij corrected time series, one for each 
value of first order spin-down. For the current analysis we 
have rij = 33 for Crab and rij — 3 for Vela, see Tab. HU) 
where also other relevant quantities are given. In partic- 
ular, the total number of points in the frequency /spin- 
down plane is about 5.28 x 10 6 for Crab analysis, and 
4.8 x 10 5 for Vela analysis. For each time series we ap- 
ply a final cleaning step by removing the largest outliers. 
These are identified by histogramming the logarithm of 


the absolute value of the data amplitude and choosing 
a threshold approximately marking the beginning of the 
non-Gaussian tail of data distribution. In Fig. [5] the his- 
togram of Crab and Vela data amplitude, correspond- 
ing to the “central” time series, are shown. For Crab a 
threshold of 1.5 x 10 — 21 , corresponding to -1.82 in the 
figure x-axis, has been used to remove outliers, while for 
Vela a value of 1.2 x 10 -20 , corresponding to -0.92, has 
been used. Correspondingly, the fraction of removed data 
is of about 1.8% for Vela and 1.9% for Crab. By apply- 
ing the Kolmogorov- Smirnov test we have verified that 
for each pulsar the data distribution of the various time 
series are fully in agreement thus justifying the use of the 
same threshold for all of them. 

At this point the detection statistic is computed for ev- 
ery frequency and spin-down value in the explored range. 
The detection statistic we use is based on the so-called 5- 
vectors, the same used for pulsar targeted searches il, 
and is here briefly described. Once barycentric and spin- 
down corrections have been applied, a CW signal with 
frequency /o present in the data would be monochro- 
matic, apart from an amplitude and phase sidereal mod- 
ulation due to the time-varying detector beam pattern 
functions, and given by Eq. m with uo(t) constant and 
equal to 27r/o- From Eqs. (I2I5I6I) it follows that the signal 
is completely described by its Fourier components at the 
5 angular frequencies cjo, <^o cjo±2f} 0 . This set of 

5 complex numbers constitutes the signal 5-vector. Given 
a generic time series g(t), the corresponding 5-vector is 
defined as 


G = J fl (t)e- t ("°- kn ®) t dt (22) 

where k = [— 2, — 1, ..., 2] and T is the observation time. 
Let us indicate with X the data 5-vector and with 
A+, A x the signal plus and cross 5- vectors, obtained by 
applying the definition of Eq. ([ 22 ]) to Eqs. (|5l6p . These 
two last quantities depend only on known parameters and 
form the signal templates. Once the 5- vectors of data and 
of signal templates have been computed, the two complex 
numbers 


X • A + / x 

^+/ x — |a+/x 1 2 

are built, see 10 for more details. They correspond to 
computing matched filters between the data and the sig- 
nal templates, and it can be shown, assuming the noise 
is Gaussian with mean value zero, that they are esti- 
mators of the signal plus and cross complex amplitudes 
Hoe l ®°H x - These estimators are used to 
build the detection statistic 

5= |A+| 4 |iJ +| 2 + |A x | 4 | J ff x | 2 (24) 

The maximum of the detection statistic over the 
searched parameter space, <S macc , is determined. This 
is the loudest candidate and is identified by a triple 
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FIG. 4. The parallelograms delimit the portion of the / — / plane covered in the narrow-band search for the Vela (upper plot) 
and Crab (bottom plot) pulsars. Only one value of / has been considered in the analysis. The total number of points is given 
in the last column of Tab. m 




FIG. 5. Histogram of the logarithm of data amplitude for Vela (top plot) and the Crab (bottom plot), used to select the threshold 
for the removal of outliers. For Crab the threshold has been put at -1.82, corresponding to an amplitude of ~ 1.5 x 10 — 21 , 
while for Vela a value of -0.92, corresponding to an amplitude of ~ 1.2 x 10 — 20 , has been chosen. 
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TABLE I. The Crab and Vela reference parameters, a and 8 are the equatorial coordinates, fo, fo, fo are the gravitational- wave 
frequency parameters. The reference epoch for the Crab position is MJD 54632, while for the rotational parameters is MJD 
55696. The reference epoch for Vela is MJD 53576 both for position and GW frequency parameters. The Crab ephemeris has 
been obtained from a fit of the Jodrell Bank monthly ephemeris [3 , while for Vela they have been derived from Hartebeesthoek 
radio telescope observations 0. Also estimations of the polarization parameters and their associated uncertainty, obtained 
from the analysis of Chandra X-ray observations of the pulsar wind nebula torus 2(1.0. are given in the last two columns. 
The analysis presented in 0,0 does not allow to determine the sense of the star’s spin, so values of i and ip corresponding 
to l — y 180° — l, ip —tip + 180° are also possible. However, the upper limits on strain amplitude reported in this paper are not 
sensitive to the sense of rotation. 


Source 

a [hh:mm:ss] 

<5 [deg] 

fo [Hz] 

fo [Hz/s] 

/o[Hz/s 2 ] 

t[deg] 

V-[deg] 

Crab 

05 : 34 : 31.97 

22.0145 

59.4448 

-7.4183 x l(T iu 

2.6307 x 10 _2U 

62°. 2 =b 1°.9 

35°. 2 =b 1°.5 

Vela 

08 : 35 : 20.61 

-45.1764 

22.3840 

-3.1460 x 10 _ii 

1.2848 x 10 _2i 

63°.6±0°.6 

40°. 6 d= 0°.l 


TABLE II. Main quantities related to the parameter space of VSR4 narrow-band search. A / is the frequency range; 8 is the 
width parameter defined in Eq. ©; A / is the first order spin-down range; n/ is the number of frequency bins; ny and rif are, 
respectively, the number of bins for the first and second order spin-down; ntot — rif - Uj is the total number of points in the 
parameter space. 


Source 

A/[Hz] 

8 

A/[Hz/s] 

U f 

n f 

n f 

Tltot 

Crab 

0.02 

1.68 x 10 -4 

2.49 x 10 -13 

1.6 x 10 s 

33 

1 

5.28 x 10 B 

Vela 

0.02 

4.47 x 10 -4 

2.81 x 10 -14 

1.6 x 10 s 

3 

1 

4.80 x 10 5 


{Smax i fSmax , fsmax ) • This means we need to compute 
ntot = nf • rij values of the detection statistic. 

The maximum value S max is used to assess detection 
significance by computing the p- value, that is the prob- 
ability that a value of the detection statistic equal, or 
larger, than Smax can be obtained in the absence of 
any signal. It implies the need to compute the noise- 
only distribution of the detection statistic. This is a 
multi-dimensional probability distribution (with dimen- 
sion n t ot ) which would be difficult to compute and to han- 
dle. In practice, as discussed in the p - value is com- 
puted considering the single-trial noise probability distri- 
bution, that is the same that would be used for a targeted 
search, and choosing a suitable threshold on it to dis- 
criminate between interesting (that is, deserving a deeper 
study) and not interesting candidates. As shown in [llj], 
by setting an overall p-value po (over the full multi- 
dimensional distribution) the corresponding threshold on 

the single-trial distribution is p t hr — 1 — (1 — po ) ntot . In 
our case, by setting, e.g., po = 0.01 the resulting p t hr 
would be of the order of 10 -8 . In principle, we would 
like to generate the noise-only probability distribution by 
computing the detection statistic at several “off-source” 
frequency bins, that is frequencies near but outside the 
explored range, where we are assuming the signal could 
be. This is what is typically done in targeted searches ( 3 ). 
In the present analysis, however, in order to appreciate p- 
values of the order of p t hr we should consider of the order 
of 10 8 off-source frequency bins. This is computationally 
impractical. For this reason we use the theoretical dis- 
tribution, which assumes the noise is Gaussian, given in 
0 : 


m 


I A X |2 


e 


r |A +| 2 


(25) 


where o\ = a 2 • T, with a 2 being the noise variance. 

If Smax is compatible with noise at a given confidence 
level, e.g. 1%, we compute an upper limit on signal 
strength, using the same method described in 4] , other- 
wise in case of detection signal parameters are estimated 
through suitable combinations of the real and imaginary 
parts of H+ and H x , as explained in 0- 


V. RESULTS 

Analysis results are summarized in Tab. IIII1 The data, 
both for the Vela and Crab searches, are compatible with 
noise. In particular, we find p - values equal to 0.33 and 
0.013 for Vela and the Crab respectively, larger than 
po = 0.01 chosen to identify interesting candidates. For 
Crab, however, the obtained p - value is very near to the 
chosen threshold. Altough a reasonable choice of po is 
rather arbitrary (we could have for instance chosen 10 -3 ), 
we have decided to study in some detail the candidate’s 
properties. In particular, we have considered the distri- 
bution of the top ten candidates in the frequency /spin- 
down plane. They appear to be randomly distributed, 
without the clustering that we would expect in presence 
of a signal. To verify this hypothesis we have added 
a simulated signal to the data, with the same parame- 
ters as the Crab, but with a slightly different frequency 
and spin-down, and with an amplitude ho ~ 4.1 x 10 -25 
such that the resulting loudest candidate has a value of 
the detection statistic approximately equal to the loud- 
est candidate of the actual analysis. We have then run 
the full analysis on this new data set and looked again at 
the top ten candidate distribution in the frequency /spin- 
down plane. In this case we indeed observe a clustering, 
with 4 out of 10 candidates having a frequency within 


°x (I A x | 2 - |A+| 2 ) 
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±2 bins of the injected value, and 4 out of 10 candidates 
having a spin-down within ±3 bins of the injected value. 
Overall, 5-6 of the top ten candidates appear to be due 
to the injected signal. Top ten candidates distribution in 
the frequency /spin-down plane is shown in Fig.[6]for both 
cases. The results of this test make us more confident in 
declaring a non-detection also in the Crab analysis. 


We have computed 95% confidence level upper limits, 
both assuming uniform priors on the polarization param- 
eters, ^ and cos^, and “restricted” priors described by 
a Gaussian distribution with mean value and standard 
deviation given in Tab. HI In fact the upper limits we 
compute are obtained from the posterior distribution of 
the signal strain amplitude, conditioned to the observed 
value of the detection statistic, as described in [4]. For 
Vela the upper limits are very similar to the spin-down 
limit. This does not allow us to significantly constrain the 
fraction of spin-down energy due to the emission of grav- 
itational waves. In the case of the Crab pulsar the upper 
limits on signal strain amplitude are about 2 times be- 
low the spin-down limit, with a corresponding constraint 
of about 25% on the fraction of spin-down energy due 
to gravitational waves. The upper limit on signal strain 
amplitude can be converted, via Eq. m, into an upper 
limit on star ellipticity of about 3.7 x 10 -4 , assuming the 
neutron star moment of inertia is equal to the canonical 
value of 10 38 kg m 2 . The upper limits on ellipticity are 
comparable to the maximum value foreseen by some “ex- 
otic” equation of state for neutron star matter [7, [1| , but 
are much larger than the maximum value predicted for 
a “standard” equations of state. In Tab. IIIII upper limits 
on the quadrupole mass moment, which are independent 
on the uncertain value of the star moment of inertia, are 
also given. An uncertainty of about 8% is associated with 
the upper limit on strain amplitude. This has been es- 
timated as the square root of the quadratic sum of the 
calibration error which, as discussed in Sec. um amounts 
to 7.5%, and of the error associated to the finite size of 
the simulation used to compute the upper limits, which 
is about 3%. 


While this is the first time a narrow-band search has 
been carried out for the Vela pulsar, our results for 
the Crab represent an improvement with respect to the 
LIGO/S5 search described in Q. Overall, this is the first 
time the spin-down limit has been significantly overcome 
in any narrow-band search. On the other hand, our up- 
per limits are clearly worse than those found in Vela and 
Crab targeted searches [4]. As explained in m this is 
what we expect as a consequence of the lower sensitivity 
of narrow-band searches, due to the volume of parameter 
space that is explored. Also, the use of VSR4 data alone 
for the Crab search contributes, even if to a lesser extent, 
to reduce the search sensitivity with respect to Q where 
data from Virgo VSR2, VSR4 and LIGO S6 were used. 


VI. VALIDATION TESTS 

We have performed validation tests of the analysis 
pipeline both injecting software simulated signals in 
Gaussian noise and performing a narrow-band search 
around some of the hardware injected signals in VSR4 
data. The first kind of test has been discussed in Q 
Here we focus attention on the second kind of test. 

For the entire duration of VSR4 run, 10 simulated CW 
signals have been injected in the Virgo detector by send- 
ing the appropriate excitations to the coils used to control 
mirror position. These signals were characterized by var- 
ious amplitudes, frequency, spin-down and polarization 
parameters, and corresponded to sources with various lo- 
cations in the sky. In particular here we have considered 
injections named PulsarS, Pulsar5, Pulsar8 with param- 
eters of signal amplitude ho, frequency /o, spin-down /o, 
position (aoj^o)? ratio between the polarization ellipse 
semi-minor and semi-major axes 770 and wave polariza- 
tion angle 7 / 0 , given in Tab.HVl For each of the hardware 
injections we have performed a narrow-band search over a 
frequency range of 10 -4 Hz and over a spin-down range of 
1.585 x 10 -13 Hz/s around the signal injected values. The 
grid in frequency and spin-down has been built in such 
a way that each center bin corresponds to a frequency 
given by an integer number of bins. This means that the 
true values do not coincide with a bin center. The ex- 
plored frequency range is covered by 814 bins while the 
spin-down range is covered by 21 bins. For each hardware 
injection the frequency and spin-down corresponding to 
the loudest candidate have been selected and compared 
to the true signal values. Moreover, the amplitude and 
polarization parameters of the injected signals have been 
estimated and compared to the actual values and to the 
values that would have been found in a targeted search. 
In Tab. El and EH test results are summarized. In par- 
ticular, in Tab. El we report the estimated frequency /, 
the error with respect to the true value e/, measured in 
number of frequency bins, the estimated spin-down value 
/ and its error ey in number of spin-down bins. From 
Tab. E| we see that, within the discretization error, both 
the frequency and the spin-down are correctly recovered 
for all the hardware injections. This can be also seen in 
Fig. □ where two plots, referring to PulsarS , are shown. 
The upper plot represents the maximum of the detection 
statistic (maximized over frequency) as a function of the 
spin-down index, which goes from 1 to 21. The maximum 
of the curve is at index number 11 which is the nearest 
to the injected value. Note that given the very small 
spin- down value of PulsarS and Pulsar 5, much smaller 
than the spin-down bin width, the searched range covers 
both negative and positive values of /. The lower plot 
shows the value of the detection statistic as a function of 
the frequency bin, which goes from 1 to 814, computed 
taking the spin-down of the loudest candidate. As be- 
fore the maximum is found in correspondence of the bin 
nearest to the injected value of the frequency. In Tab. I VII 
the error in the estimation of amplitude and polarization 
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FIG. 6. Distribution in the frequency /spin- down plane of the top ten candidates obtained in the actual analysis (blue circles) 
and obtained after the addition to the data of a simulated signal with amplitude ho = 4.06 x 10 -25 , such that the resulting 
loudest candidate is approximately as loud as the one obtained in the actual analysis (red squares) . The black star identify the 
injection. Frequency and spin-down are expressed as number of bins from the beginning of the corresponding intervals. Six red 
squares appear to be at the same frequency, but this is just a visual effect due to the large range of frequency covered by the 
x-axis. In fact, they are very near but spread around the injection values. 


TABLE III. Summary of the analysis results. The p-value is reported in the second column; h^J 1 ^ is the experimental upper 
limit assuming uniform priors for the polarization parameters; h^f[ tr is the experimental upper limit assuming “restricted” 
priors for the polarization parameters, given by a Gaussian distribution with mean value and standard deviation given in 
Tab. [2 h s d is the spin-down limit, computed through Eq. m; euL is the upper limit on the ellipticity; Q 22 ,ul is the upper 
limit on the mass quadrupole moment; huL/h s d is the ratio between the upper limit and the spin-down limit on signal strain 
amplitude; EuL/E s d is the upper limit on the fraction of spin-down energy due to the emission of gravitational waves (the 
values among parentheses refer to the restricted case). 


Source 

p- value 

lunif 

,l UL 

l restr 
n UL 

hsd 

cul x 10 4 

Q 22 ,ul X 10 34 [kg m 2 ] 

huL/h s d 

EuL/E s d 

Vela 

0.33 

3.2 x IQ-' M 

3.3 x 10~' Z4 

3.3 x 10~ 24 

17.6(19.1) 

13.4(13.9) 

0.97(1) 

0.94(1) 

Crab 

0.013 

7.0 x 10“^ 

6.9 x 10“^ 

1.4 x W™ 

3. 8(3. 7) 

2. 9(2. 8) 

0.50(0.49) 

0.25(0.24) 


parameters of hardware injections is given. In particular, 
in column 2 the ratio between the found and the injected 
amplitude is reported, in column 3 the fractional error 
on 77 , that is e v = and in column 4 the fractional 

error on t/j, . Values in parentheses refer to a 

targeted search for the same injections. Parameters are 
generally well recovered, with an accuracy just slightly 
worse with respect to the targeted search case. This is 
due to the small error in the estimation of the signal fre- 
quency in a narrow-band search, consequence of the finite 
size of the grid step. In fact, looking for instance at the 
amplitude estimation, we can note that for Pulsar8 , for 


which the frequency error is the largest, the loss with re- 
spect to the targeted search is the biggest, as expected. 


VII. CONCLUSIONS 

Targeted searches for continuous gravitational wave 
signals assume a strict correlation between the gravita- 
tional wave frequency and the star rotation frequency. 
For instance, for a neutron star non-axisymmetric with 
respect to the rotation axis the gravitational wave fre- 
quency is exactly two times the rotation rate. However, 
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TABLE IV. Main parameters of the hardware injections used to test the narrow-band search pipeline, ho is the signal amplitude; 
fo is the signal frequency; fo is the spin-down, (ao, £o) is the source position in equatorial coordinates, 770 is the ratio between 
the polarization ellipse semi-minor and semi-major axes, ipo is the wave polarization angle. 


Name 

h 0 

fo [Hz] 

/ o[Hz/s] 

«o[deg] 

<5o[deg] 

Vo 

i>o [deg] 

Pulsar 3 

8.296 x 10~ 24 

108.857159396 

-1.46 x 10“ iv 

178.372574 

-33.436602 

0.1600 

25.439 

Pulsar 5 

3.703 x 10 “^ 4 

52.808324359 

-4.03 x 10 " 18 

302.626641 

-83.8391399 

-0.7625 

-20.853 

Pulsar 8 

8.067 x 10 - * 4 

192.237058812 

-8.65 x 10 “ 9 

351.389582 

-33.4185168 

-0.1470 

9.7673 


TABLE V. Tests with hardware injections: frequency and spin-down recovery. For each hardware injection the estimated 
frequency and spin-down, / and /, and the associated errors e/ and ey, expressed in number of bins, are given. 


Name 

/[Hz] 

e f 

/[Hz/s] 

6 f 

Pulsar 3 

108.857159411 

0.12 

0 

0.002 

Pulsar5 

52.808324371 

0.10 

0 

0.0005 

Pulsar 8 

192.237058796 

-0.13 

-8.6500021 x 10 -9 

-0.272 



FIG. 7. Tests with hardware injection PulsarS. Upper plot: maximum of the detection statistic (maximized over frequency) as 
a function of the spin-down index (which goes from 1 to 21). Bottom plot: value of the detection statistic as a function of the 
frequency bin (which goes from 1 to 814), computed taking the spin-down of the loudest candidate. The plots show that the 
frequency and spin-down of the loudest candidate are the nearest to the true values. 


TABLE VI. Tests with hardware injections: amplitude and polarization parameter recovery. For each hardware injection we 
report the ratio between the estimated and injected signal amplitude, h/ho, the relative error on 77 normalized to its range of 
variation (2), e^, and the relative error on 7/; normalized to its range of variation (90 degrees), e^. The numbers in parentheses 
refer to a targeted search for the same signals. 


Name 

h/ho 


! 1 

Pulsar3 

1.016(0.992) 

0.0027(0.0021) 

1 

O 

O 

O 

T 

0 

0 

0 

OO 

Pulsar5 

1.023(1.002) 

0.0080(0.0089) 

0.1196(0.0104) 

Pulsar 8 

0.871(0.975) 

0.0278(0.0038) 

0.0007(0.0064) 
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it is questionable that such strict a correlation is always 
valid or that is maintained over long times, and various 
mechanisms have been proposed that could produce a 
mismatch. For this reason it is important to have in place 
an analysis procedure robust with respect to deviations 
from the standard assumption made in targeted searches 
and then able to perform a narrow-band search for CW 
signals over a range of frequency and spin-down values. 
In this paper we have presented results of a narrow-band 
search for CW from the Crab and Vela pulsars in the 
data of Virgo VSR4 run. In both cases the data appear 
to be fully compatible with noise, and we have set 95% 
confidence level upper limits on strain amplitude assum- 
ing both uniform and restricted priors on polarization 
parameters. The upper limits are, respectively, compa- 
rable (for Vela) and a factor of 2 below (for the Crab) 
to the corresponding spin-down limits and while this is 
the first time a narrow-band search has been done for the 
Vela pulsar, for the Crab our results significantly improve 
with respect to past analyses. As expected, the narrow- 
band search upper limits are worse than those established 
in targeted searches of the same sources. This analysis 
method will be applied to narrow-band searches of CW 
signals from several potentially interesting pulsars in data 
of advanced Virgo and LIGO detectors, which will start 
to collect science data at the end of 2015, and are ex- 
pect to reach their target sensitivity, about one order of 
magnitude better than first generation detectors, around 
2018 
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